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Abstract 

Explaining cases of long-term persistence of parthenogenesis has proven an arduous task for evolutionary biologists. Interpreting sex¬ 
ual-asexual interactions though has recently advanced owing to methodological design, increased taxon sampling and choice of model 
organisms. We inferred the phylogeny of Artemia, a halophilic branchiopod genus of sexual and parthenogenetic forms with cosmopol¬ 
itan distribution, marked geographic patterns and ecological partitioning. Joint analysis of newly derived ITS1 sequences and 16S RFLP 
markers from global isolates indicates significant interspecific divergence as well as pronounced diversity for parthenogens, matching that 
of sexual ancestors. Maximum parsimony, maximum likelihood, and Bayesian methods were largely congruent in reconstructing the phy¬ 
logeny of the genus. Given the current sampling, at least four independent origins of parthenogenesis are deduced. Molecular clock cal¬ 
ibrations based on biogeographic landmarks indicate that the lineage leading to A. persimilis diverged from the common ancestor of all 
Artemia species between 80 and 90 MYA at the time of separation of Africa from South America, whereas parthenogenesis first appeared 
at least 3 MYA. Common mitochondrial DNA haplotypes delineate A. urmiana and A. tibetiana as possible maternal parents of several 
clonal lineages. A novel topological placement of A. franciscana as a sister clade to all Asian Artemia and parthenogenetic forms is pro¬ 
posed and also supported by ITS1 length and other existing data. 

© 2006 Elsevier Inc. All rights reserved. 
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1. Introduction 

Evolution has solved the problem of reproduction in 
various ways. Strictly asexual organisms though—those 
that reproduce without fertilization—are, in the long-term, 
something of an evolutionary paradox. They are barred 
not only from the benefits of sexual reproduction (out- 
crossing and meiotic recombination) but also from the 
resulting genetic variance needed for evolutionary change. 
Current theories attempting to explain the dominance of 
sex (Judson and Normark, 1996; Normark et ah, 2003) 
are largely unsuccessful in accounting for parthenogenetic 
lineages. Most notably, the inability of such taxa to purge 
deleterious mutations is thought of as the prime factor 
driving them to evolutionary degradation (Hurst and Peck, 
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1996). Nevertheless, parthenogens attain short-term persis¬ 
tence and, occasionally, provocative longevity (Schon 
et al., 1996; Mark Welch and Meselson, 2000). For this rea¬ 
son, the origin, genetic diversity and phylogenetic history 
of asexual systems have been atop evolutionary biologists’ 
list of inquiries. 

Mutational and ecological models provide more or less 
specific predictions regarding the genetic and life history 
architecture of unisexual taxa, the patterns of variation in 
nuclear and organelle DNA and the timing of consequenc¬ 
es due to loss of sexuality (Normark and Moran, 2000; 
Normark et al., 2003; Simon et al., 2003). In addition, 
the rate and mode of origin of parthenogenesis may be cru¬ 
cial in determining levels of genetic diversity, the geograph¬ 
ic distribution and ecological success of asexuals (Sandoval 
et al., 1998; Vrijenhoek, 1998; Schon et al., 2000; Law and 
Crespi, 2002; Paland et al., 2005). Phylogenetic approaches 
to the evolution and lifespan of asexuality are especially 
powerful since typical signatures of the loss of sex 
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(Delmotte et al., 2001), rates of molecular change (Schon 
et al., 2003) and hybridization (Delmotte et al., 2003) can 
be readily detected on reconstructed trees. Methodological¬ 
ly, such investigations are largely dependent on adequate 
sampling of extant taxa and sufficient screening of genomic 
regions. The former may seriously influence both the esti¬ 
mated timeframe of asexual lineages and the number of 
transition events. The latter can help determine whether 
the observed clonal divergence is mutationally generated 
or recurrently captured from sexual ancestors (Crease 
et al., 1989; Chaplin and Hebert, 1997). It is thus obvious 
that more robust interpretations are made from model sys¬ 
tems where the alternative modes of reproduction are con¬ 
founded within a single organismal lineage and additional 
biogeographic and ecological data can be related to pat¬ 
terns of divergence. Recent work has emphasized the need 
for this integration (Simon et al., 2003). It is paradoxical 
though that, for Artemici , admittedly the most accom¬ 
plished survivor of hypersaline settings and a confusing 
case of ancient asexuality, similar critical and detailed 
assessments are overdue. 

Artemia is a genus (Crustacea, Anostraca) of sexual and 
parthenogenetic forms with a global distribution in inland 
salt lakes and coastal lagoons (Triantaphyllidis et al., 1998; 
Van Stappen, 2002). The organism has featured in the liter¬ 
ature extensively, by virtue of its importance in aquaculture 
and as a model system for varied research. From an evolu¬ 
tionary perspective in particular, different adaptations, dis¬ 
tinctive genetic features and marked biogeographic 
patterns are all found in the genus, thus offering unique 


opportunities for studies on phylogeny and the interaction 
between sexuality and parthenogenesis (Abatzopoulos 
et al., 2002a). Currently, six bisexual species (Abatzopoulos 
et al., 2002b) and a heterogeneous group of obligate par- 
thenogens are recognized, which are either apomictic or 
automictic. A sharp geographic boundary separates the 
New World bisexuals (A. franciscana , A. persimilis ) from 
their Old World relatives (A. salina, A. urmiana , A. sinica, 
and A. tibetiana). Similarly, parthenogenetic populations 
are restricted to the Old World, where they comprise the 
majority. Bisexuals are diploid with 2 n = 42 (except A. per¬ 
similis where 2 n = 44), while parthenogens range in ploidy 
from 2 n to 5 n (Abatzopoulos et al., 1986). Morphometric 
and/or morphological, life history and genetic divergences 
are widely partitioned both within and between the differ¬ 
ent reproductive modes (for review see Browne, 1992). 
Allopatric divergence and ecological specialization are 
believed to have shaped Artemia evolution, while the influ¬ 
ence of dispersal on contemporary regional distributions 
has only recently been explored (Green et al., 2005). 

A historical aspect of Artemia phylogeny has been main¬ 
ly obtained through allozyme studies. Beardmore and 
Abreu-Grobois (1983) outlined the series of phylogenetic 
events in the genus. Evaluation of the degree and patterns 
of interspecific divergence have indicated that the primal 
evolutionary event has been the separation of New and 
Old World bisexual lineages. This was followed by the sep¬ 
aration of A. franciscana and A. persimilis in the New 
World and the divergence of A. salina and A. urmiana lines 
in the Old World (see tree in Fig. 1). Based on an allozymi- 



Fig. 1. The 46 Artemia individuals (sequenced for ITS1 region) as distributed in the world map. O, A. franciscana ; O. A. persimilis ; □, A. salina ; A, A. 
urmiana ; +, A. tibetiana ; -f. A. sinica ; and A parthenogenetic Artemia strains. Inserted tree shows currently accepted phylogenetic relationships in the 
genus as determined by allozyme markers (see also Beardmore and Abreu-Grobois, 1983). Numbers at nodes are divergence times in million years. 
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cally calibrated molecular clock, Abreu-Grobois (1987) 
speculated on a Mediterranean genesis of parthenogenetic 
Artemia. Accordingly, the parthenogenetic clade branched 
from an Old World sexual ancestor (possibly A. urmiana) 
approximately 5.5 million years ago (MYA). Unreliable 
as it may be from a molecular marker perspective (see 
Avise, 1994) this estimate strikingly accords with geological 
evidence suggesting a synchronous (5.96-5.33 MYA, 
Krijgsman et ah, 1999) dramatic increase in salinity and 
subdivision of habitats in the Mediterranean basin (Hsu, 
1972). In those conditions, selection favoured enhanced 
colonization potential granted by reproductive assurance 
(i.e. parthenogenetic reproduction, Baker, 1965). Subse¬ 
quent events involved rapid radiation of asexual forms in 
a wide range of habitats coupled with ploidy shifts. More 
recently though, a substantially more ancient origin of par¬ 
thenogenetic Artemia (30—40 MYA) was claimed by Perez 
et al. (1994) based on mtDNA sequence divergence. 

Despite extensive recording of Artemia biodiversity 
(Abatzopoulos et ah, 2002a) and recurrent evidence for 
eccentric sequence change (Maley and Marshall, 1998; 
Crease, 1999; Saez et ah, 2000; Coleman et ah, 2001; 
Hebert et ah, 2002), full-scale phylogenetic analyses of 
the genus are lacking. Based on a global sampling and 
using sequences of the ITS1 nuclear region, we have 
tried to infer the phylogenetic relationships and investi¬ 
gate the biogeographic history and ages of all recognized 
Artemia species and asexual forms. To assess congruence 
between nuclear and mtDNA markers, a subset of pop¬ 
ulations were also screened for 16S mtDNA restriction 
patterns. 

Throughout this study we have chosen to refer to asex¬ 
ual forms of Artemia as populations. The indiscreetly used 
term A. parthenogenetica is at odds with the biological spe¬ 
cies concept, since obligate parthenogens existing in numer¬ 
ous clones at several different ploidy levels cannot readily 
be considered as belonging to a single species (Abatzopou¬ 
los et ah, 2002a). 

2. Materials and methods 

2.1. Populations analyzed 

Artemia exhibits two types of parturition. Embryos in 
the brood pouch of the female either develop to free-swim¬ 
ming nauplii (ovoviviparity) or are released as cysts (ovi- 
parity), forming resting propagule banks which may 
sequester a substantial proportion of variability from the 
active population. Therefore, due to the significance of 
resting cyst banks to evolutionary change (Hairston and 
De Stasio, 1988; De Meester et ah, 2002) and to sample 
overall genetic diversity more adequately, we have used 
encysted embryos as starting material for the analyses. 
Sampling was on a global scale (Fig. 1) and included all 
major biogeographical regions and, when possible, multi¬ 
ple representatives from the different localities (Table 1). 
Cyst hatching (hatching percentage >85%) and rearing of 


nauplii (survival >90%) followed procedures described in 
Baxevanis et ah (2005). Samples with lower values for these 
two indices were excluded. With the onset of sexual matu¬ 
rity, examination of sex ratios enabled the unambiguous 
confirmation of parthenogens, in conjunction with pub¬ 
lished data on these populations. 

2.2. DNA extraction, PCR, and sequencing 

DNA was extracted from single specimens (starved for 
at least 24 h) using 150-200 pi of Chelex suspension (6%, 
Bio-Rad Laboratories, CA, USA) with proteinase K 
(15 pi of 20mg/ml initial stock) pre-treatment (Estoup 
et ah, 1996). “No-DNA” extractions were also conducted 
to check for contamination. PCR and sequencing of the 
complete ITS1 region was performed using the universal 
primer 18d-5' by Hillis and Dixon (1991) and a newly 
designed primer R58 (5' ACCC ACGAGCC AAGT G 
ATCC 3') based on deposited 5.8S rDNA sequences of 
Penaeus vannamei (AF121132). Amplifications were per¬ 
formed on a PTC-100® Peltier thermal cycler (MJ 
Research) with the following temperature profiles and con¬ 
ditions: 4 min at 93 °C, 32 cycles of 40 s at 93 °C, 40 s at 
62 °C, 1 min at 72 °C and a final extension of 5 min at 
72 °C. Total reaction volumes of 25 pi consisted of 2.5 pi 
template DNA, 5 pi 5x PCR buffer, 2.5 mM MgCl 2 , 
1 mM dNTPs, 0.001 mM of each primer and 1.25 U of 
Taq DNA polymerase (Expand High Fidelity PLUS PCR 
System, Roche). PCR products were excised from agarose 
gels (2%) and purified using the NucleoSpin Extract kit 
(Macherey-Nagel). Cycle sequencing was performed with 
Excell II DNA sequencing kit (EPICENTRE) in both 
directions. Due to the size of the PCR product (see Section 
3), an additional set of internal primers were also designed 
(sequences available upon request) and used for sequenc¬ 
ing. Reactions were electrophoresed on a PRISM 3730x1 
DNA analyzer (Applied Biosystems). 

2.3. Sequence alignment and phylogenetic analyses 

The total data set consisted of 47 sequences (outgroup 
included) (see Table 1). BLAST searches confirmed the 
identity of the PCR product on the basis of detected simi¬ 
larities with 18S and 5.8S rRNA genes from a number of 
taxa as well as A. franciscana AJ238061, A. salina 
X01723 (see Weekers et ah, 2002), Artemia sp. M33097 
(see Vaughn et ah, 1984). A fairy shrimp sequence ( Strep- 
tocephalus proboscideus, AY519840, 680 bp), spanning part 
of 18S, the complete internal transcribed spacer 1 and part 
of 5.8S (Daniels et ah, 2004), was used to root the phylog- 
eny (used sequence length = 612 bp). For alignment we 
used the program MEME (Multiple Em for Motif Elicita¬ 
tion, Bailey and Elkan, 1994) in combination with 
MACAW (Schuler et ah, 1991) to identify regions of statis¬ 
tically significant local similarity. In total, six highly con¬ 
served regions in both the ingroup taxa and the outgroup 
were chosen to guide the final alignment. Two of them were 
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Table 1 

List of the studied individuals, with their abbreviations, species designation and GenBank accession numbers for ITS1 sequences 


Individual 

Population 

Species 

Accession No. 

SFB1 

San Francisco Bay, USA 

A. franciscana a 

DQ201298 

SFB2 

San Francisco Bay, USA 

A. franciscana a 

DQ201297 

SFB5 

San Francisco Bay, USA 

A. franciscana a 

DQ201293 

SFB8 

San Francisco Bay, USA 

A. franciscana a 

DQ201292 

SFB9 

San Francisco Bay, USA 

A. franciscana a 

DQ201294 

GSL1 

Great Salt Lake, USA 

A. franciscana 3 

DQ201300 

GSL3 

Great Salt Lake, USA 

A. franciscana 0 

DQ201301 

IQU1 

Iquique, Chile 

A. franciscana L 

DQ201296 

ING3 

Ingebright, Canada 

Unknown 15 

DQ201299 

MAS1 

Margherita di Savoia, Italy 

A. franciscana d 

DQ201295 

MAS6 

Margherita di Savoia, Italy 

A. franciscana ' d 

DQ201291 

MAS7 

Margherita di Savoia, Italy 

A. franciscana d 

DQ201290 

BAI7 

Buenos Aires, Argentina 

A. persimilis a 

DQ201263 

BAI8 

Buenos Aires, Argentina 

A. persimilis a 

DQ201265 

BAI9 

Buenos Aires, Argentina 

A. persimilis a 

DQ201264 

PIC7 

Pichilemu, Chile 

A. persimilis c 

DQ201266 

TDP2 

Torres del Paine, Chile 

A. persimilis c 

DQ201268 

TDP3 

Torres del Paine, Chile 

A. persimilis c 

DQ201267 

SFA4 

Sfax, Tunisia 

A. salina a 

DQ201302 

SFA9 

Sfax, Tunisia 

A. salind a 

DQ201305 

SAN4 

Sant’ Antioco, Italy 

A. salinct 

DQ201304 

SAN6 

Sant’ Antioco, Italy 

A. salina e 

DQ201303 

WNA4 

Wadi El-Natrun, Egypt 

A. salincf 

DQ201306 

WNA5 

Wadi El-Natrun, Egypt 

A. salina 

DQ201307 

LAR8 

Larnaca, Cyprus 

A. salincr 

DQ201308 

LAR9 

Larnaca, Cyprus 

A. salina 3 

DQ201309 

URM1 

Urmia Lake, Iran 

A. urmiana a 

DQ201276 

URM2 

Urmia Lake, Iran 

A. urmiana a 

DQ201277 

URM3 

Urmia Lake, Iran 

A. urmiana a 

DQ201275 

SIN1 

Yuncheng, P.R. China 

A. sinica a 

DQ201285 

SIN2 

Yuncheng, P.R. China 

A. sinica a 

DQ201286 

SIN3 

Yuncheng, P.R. China 

A. sinica a 

DQ201287 

TIB3 

Lagkor Co, Tibet, P.R. China 

A. tibetiana 

DQ201270 

TIB7 

Lagkor Co, Tibet, P.R. China 

A. tibetiana i 

DQ201269 

MEM4 

Megalon Embolon, Greece 

Parthenogenetic b 

DQ201274 

MEM9 

Megalon Embolon, Greece 

Parthenogenetic b 

DQ201273 

JIA7 

Jiangsu, P.R. China 

Parthenogenetic b 

DQ201288 

POL14 

Polychnitos, Greece 

Parthenogenetic 8 

DQ201283 

EIL1 

Eilat, Israel 

Parthenogenetic b 

DQ201272 

EIL2 

Eilat, Israel 

Parthenogenetic* 5 

DQ201271 

TCL1 

Torre Colimena, Italy 

Parthenogenetic* 1 

DQ201279 

TCL2 

Torre Colimena, Italy 

Parthenogenetic* 1 

DQ201278 

MAD4 

Ankiembe, Madagascar 

Parthenogenetic 1 

DQ201280 

MAD6 

Ankiembe, Madagascar 

Parthenogenetic 1 

DQ201284 

NAM1 

Swakopmund, Namibia 

Parthenogenetic 1 

DQ201281 

NAM5 

Swakopmund, Namibia 

Parthenogenetic 1 

DQ201282 

OUTGROUP 


Streptocephalus proboscideuA 

AY519840 


a Baxevanis et al. (2005). 
b Van Stappen (2002) and references therein. 
0 Gajardo et al. (2004). 
d Mura et al. (2004). 
e Mura and Brecciaroli (2004). 

1 Abatzopoulos et al. (2002b). 

8 Triantaphyllidis et al. (1993). 
h Mura et al. (2005). 

1 Triantaphyllidis et al. (1996). 

* Daniels et al. (2004). 


located within the partial regions of 18S and 5.8S rRNA 
genes. The obtained similarities in these two regions are 
justified by the fact that the primers used in our study 
and those of Daniels et al. (2004) (for the outgroup 


sequence) amplify a similar region (partial 18S, ITS1, par¬ 
tial 5.8S). Subsequently, sequences were loaded into Clu- 
stalX 1.8 (Thompson et al., 1997) and several alignments 
were performed under a range of gap opening (7-20) and 





728 


A.D. Baxevanis et al. / Molecular Phylogenetics and Evolution 40 (2006) 724-738 


extension (1.5-12) parameters aiming to find an optimized 
multiple alignment, which preserved the structure and posi¬ 
tion of the previously identified conserved regions. This 
yielded a final alignment in ClustalX (slow-accurate mode) 
with a gap opening penalty of 15, a gap extension penalty 
of 2 and a delay divergent sequence percentage of 15. All 
sequences from this study were deposited in GenBank 
(see Table 1). 

Phylogenetic analysis was implemented in PAUP* 
4.0bl0 (Swofford, 1998) under maximum parsimony 
(MP) and maximum likelihood (ML) and in MrBayes 3.1 
(Huelsenbeck and Ronquist, 2001) under Bayesian meth¬ 
ods. We have tried not to dismiss from our analysis the his¬ 
torical information contained within gaps (see Giribet and 
Wheeler, 1999). For cladistic analysis, gaps were treated 
both as missing data and as a fifth state (according to Bena 
et al., 1998). In a third approach, they were also coded 
using GapCoder (Young and Flealy, 2003) according to 
the proposed method (“simple indel coding”) by Simmons 
and Ochoterena (2000) under the assumption that adjacent 
gap positions are caused by a single indel event. In this con¬ 
text, each indel of a particular length and position was cod¬ 
ed as a separate binary character (presence/absence) thus 
creating a matrix (266 coded indels in total) which was 
implemented along the original data matrix into phyloge¬ 
netic analyses. For MP, trees were generated using heuristic 
searches with TBR (tree-bisection-reconnection) branch 
swapping and 500 random taxon additions. Nodal support 
was assessed by 1000 bootstrap replicates. For ML (heuris¬ 
tic search; TBR branch swapping; and 100 bootstrap repli¬ 
cates) the best-fit substitution model employed was the 
HKY85 + G (Hasegawa et al., 1985), as determined by 
Modeltest 3.7 (Posada and Crandall, 1998). The parame¬ 
ters of this model were: unequal base frequencies 
(A = 0.24, C = 0.20, G = 0.28, T = 0.28), number of sub¬ 
stitution types Nst = 2, T s /T v = 0.99, y-shape parameter 
a =1.94. For Bayesian analysis (two simultaneous runs; 
1,000,000 generations; four chains; sampling every 100 gen¬ 
erations; first 2000 trees discarded as burn-in) the 
HKY85 + G model was implemented by adjusting likeli¬ 
hood settings to Nst = 2 and rates = y. For model-based 
phylogenetic methods gaps were treated as missing data. 

The assumption of molecular rate constancy was tested 
according to Takezaki et al. (1995). To infer dates of split¬ 
ting events, molecular clock calibrations were based on 
well-supported major as well as regional biogeographic 
landmarks. In reconstructing the history of the loss of 
sex, we have assumed that sexual reproduction is the ances¬ 
tral state and transitions to asexuality are irreversible. 

2.4. RFLP screening 

Total genomic DNA was separately extracted from 30 
individuals (for each population, individuals sequenced 
for ITS1 were also scored for mtDNA analysis) from the 
populations SFB, BAI, SFA, SAN, WNA, LAR, URM, 
SIN, TIB, MEM, JIA, POL, EIL, TCL, MAD, and 


NAM (abbreviations of populations are given in Table 1) 
using the procedures described previously. Part of the 
16S rRNA gene was amplified with the universal primers 
L 2 5 io and H 30 8o (Palumbi, 1996). The protocols for PCR, 
including reaction conditions and amplification profiles, 
are described in Baxevanis and Abatzopoulos (2004). 
PCR products were digested with nine restriction endonu¬ 
cleases (Alul, Bfal , Ddel, Dpnll, Haelll , Mspl, Notl, Rsal , 
and Tcuf'Y), electrophoretically separated in 1.5% agarose 
gels, stained with ethidium bromide and photographed 
under UV light. Data were coded as composite mtDNA 
genotypes based on site differences, inference of which 
was aided by the published sequence of the whole mtDNA 
of A. franciscana (Valverde et al., 1994). Nucleotide diver¬ 
gence between populations and their component haplo- 
types was computed in REAP 4.0 (McElroy et al., 1992) 
according to Nei and Tajima (1981), Nei (1987) and Nei 
and Miller (1990). A Neighbor Joining phylogram was con¬ 
structed using the Neighbor module of the PHYLIP 3.6b 
package (Felsenstein, 2004). 

3. Results 

3.1. ITS1 sequence characteristics 

The ITS1 data set for the ingroup sequences ranged 
from 983 to 1438 bp (excluding gaps in the alignment). 
The alignment length (including the outgroup) was 
1579 bp. Transition/transversion versus divergence plots 
(not shown) did not show any signs of sequence saturation. 
The disparity index statistic (Kumar and Gadagkar, 2001) 
was used to test for similarity in nucleotide composition 
across sequences and identify rogue taxa with atypical pat¬ 
terns of substitution. Such patterns were most conspicuous 
for WNA5, EIL1, EIL2, ING3, and SIN3 sequences. 
Regarding the length of the whole ITS1 region, a clade-spe- 
cific pattern was revealed with three well-defined groups. 
The first corresponded to A. persimilis sequences (bearing 
a length of 983-984 bp), the second to A. salina (1436 
1438 bp) while the third was common for A. franciscana, 
A. sinica, A. urmiana, A. tibetiana and all parthenogens 
(1133-1181 bp). This characteristic length variation was 
also further verified through sequencing of a number of 
representatives from each taxon with a second set of inter¬ 
nal primers. 

3.2. Nuclear sequence variation 

Uncorrected pairwise distances ranged from zero to 
28.1% for the ingroup taxa. Minimum and maximum 
uncorrected pairwise divergence values for the outgroup 
were 39.9 and 44.3%, respectively. Using the HKY85 + G 
model, the same estimates reached a maximum of 39.5% 
with an overall mean of 15.2% (83.8 and 17.5% including 
the outgroup, respectively). A summary of mean inter- 
and intraspecific sequence divergences is shown in Table 
2. The interspecific maximum is 34.7% for A. persimilis vs 
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Table 2 

Average inter- and intraspecific ITS1 sequence divergence estimates based on the HKY85 + G model of nucleotide substitution 



A. persimilis (3) 

A. franciscana (5) 

A. salina (4) 

A. sinica (1) 

A. urmiana (1) 

A. tibetiana (1) 

Parthenogens (7) 

A. persimilis (3) 

— 

0.242 

0.254 

0.227 

0.215 

0.217 

0.224 

A. franciscana (5) 

0.327 

— 

0.163 

0.113 

0.111 

0.108 

0.115 

A. salina (4) 

0.347 

0.199 

— 

0.132 

0.135 

0.129 

0.137 

A. sinica (1) 

0.300 

0.128 

0.154 

— 

0.031 

0.026 

0.031 

A. urmiana (1) 

0.280 

0.126 

0.159 

0.032 

— 

0.008 

0.016 

A. tibetiana (1) 

0.282 

0.123 

0.150 

0.026 

0.008 

— 

0.012 

Parthenogens (7) 

0.293 

0.132 

0.162 

0.032 

0.017 

0.012 

— 

Intraspecific divergence 

0.013 

0.018 

0.024 

0.004 

0.006 

0.000 

0.020 


Uncorrected estimates are shown above ihe major diagonal. The number of populations examined per species is indicated in parentheses. 


A. salina while the minimum is 0.8% for A. tibetiana vs A. 
urmiana. Within species divergence estimates ranged from 
zero (A. tibetiana ) to 2.4% (A. salina). Sequence diversity 
within parthenogenetic Artemia (2.0%) is considerably 
higher compared with all remaining bisexuals (except A. 
salina). In fact, certain isolates from Israel (E1L1 and 
EIL2) and China (JIA7) display divergence values from 
other asexuals (3.6-7.2%, data not shown) that greatly 
exceed interspecific comparisons among A. sinica, A. urmi¬ 
ana and A. tibetiana (see Table 2). 

3.3. MtDNA divergence 

The PCR-amplified 16S rRNA fragment was 536 bp. 
The rntDNA analysis only included populations with cor¬ 
responding ITS1 sequences. RFLP screening revealed a 
total of 30 cleavage sites, corresponding to 122 nucleotides 
(restriction patterns are available upon request). Overall, 
nine composite genotypes were identified. Sequence diver¬ 
gence estimates between populations and between haplo- 
types (designated as A-I) are given in Table 3. 
Nucleotide divergence estimates between populations ran¬ 


ged from zero to 20.2%. Similar to the results of the nuclear 
assay, A. persimilis is the most divergent lineage separated 
from all others by a distance range between 13 and 20.2%. 
Artemia urmiana and A. tibetiana are mitochondrially 
indistinguishable, sharing the same and most abundant 
haplotype (C) with five parthenogenetic populations 
(EIL, TCL, MAD, NAM, and POL). Parthenogens from 
POL differed from the previous cluster by a mean diver¬ 
gence of 1%, owing to the presence of a second haplotype 
(interclonal distance for C vs G is 5.6%, see Table 3). Par¬ 
thenogens from MEM also displayed two haplotypes (D 
and F), one of them shared with A. sinica (F). 

3.4. Phylogenetic inference 

The full 1TS1 data set (including the outgroup) consisted 
of 448 parsimony informative sites when gaps were treated 
as missing data (544 with indel coding). Reconstructed 
trees using MP, ML, and Bayesian methods were generally 
congruent and variably resolved for shallow branches. The 
50% majority rule MP tree based on indel coding (number 
of steps 1093) is shown in Fig. 2a, while that based on gaps 


Table 3 

Mean nucleotide divergence values among populations (sample size = 30 individuals) and their component haplotypes (A-I) based on 16S RFLP profiles 
with nine restriction endonucleases 


Populations 



mem d - f 

jia e 

TCL c 

pol cg 

LAR h 

EIL C 

NAM C 

MAD C 

WNA h 

SFB b 

BAI a 

SIN f 

URM C 

TIB C 

SFA 1 

SAN ] 








Composite genotypes 








mem df 

— 
















jia e 

0.037 

— 





I 

H 

G 

F 

E 

D 

C 

B 

A 


TCL c 

0.022 

0.051 

— 




— 

0.009 

0.112 

0.121 

0.131 

0.121 

0.111 

0.101 

0.174 

I 

POL cg 

0.007 

0.024 

0.010 

— 




— 

0.101 

0.117 

0.121 

0.117 

0.101 

0.092 

0.202 

H 

LAR 11 

0.103 

0.121 

0.101 

0.087 

— 




— 

0.011 

0.022 

0.035 

0.056 

0.072 

0.147 

G 

EIL C 

0.022 

0.051 

0.000 

0.010 

0.101 

— 




— 

0.010 

0.046 

0.066 

0.082 

0.157 

F 

NAM C 

0.022 

0.051 

0.000 

0.010 

0.101 

0.000 

— 




— 

0.056 

0.051 

0.092 

0.134 

E 

MAD C 

0.022 

0.051 

0.000 

0.010 

0.101 

0.000 

0.000 

— 




— 

0.020 

0.082 

0.157 

D 

WNA" 

0.103 

0.121 

0.101 

0.087 

0.000 

0.101 

0.101 

0.101 

— 




— 

0.101 

0.142 

C 

SFB b 

0.074 

0.092 

0.101 

0.074 

0.092 

0.101 

0.101 

0.101 

0.092 

— 




— 

0.101 

B 

BAI a 

0.148 

0.134 

0.142 

0.130 

0.202 

0.142 

0.142 

0.142 

0.202 

0.165 

— 




— 

A 

sin f 

0.027 

0.010 

0.066 

0.028 

0.112 

0.066 

0.066 

0.066 

0.112 

0.082 

0.157 

— 





URM C 

0.022 

0.051 

0.000 

0.010 

0.101 

0.000 

0.000 

0.000 

0.101 

0.101 

0.142 

0.066 

— 




TIB C 

0.022 

0.051 

0.000 

0.010 

0.101 

0.000 

0.000 

0.000 

0.101 

0.101 

0.142 

0.066 

0.000 

— 



SFA 1 

0.113 

0.131 

0.110 

0.097 

0.009 

0.110 

0.110 

0.110 

0.009 

0.101 

0.174 

0.121 

0.110 

0.110 

— 


SAN 1 

0.113 

0.131 

0.110 

0.097 

0.009 

0.110 

0.110 

0.110 

0.009 

0.101 

0.174 

0.121 

0.110 

0.110 

0.000 

— 
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treated as missing data (number of steps 930) is shown in 
Fig. 2b. Artemia persimilis representatives cluster into 
two groups, with sequences TDP2 and TDP3, from the 
Chilean side of the Andes, being more separated from the 
remaining Argentinean populations (bootstrap value 
100%). Clear separation of A. salina is evident into eastern 
(WNA4, WNA5, LAR8, and LAR9) and western (SFA4, 
SFA9, SAN4, and SAN6) populations with a nodal sup¬ 
port of 100%. Indel coding (Fig. 2a) has yielded a more 
robust separation of the A. salina clade (bootstrap value 
100%) as opposed to the MP tree based on missing data 
(84%, Fig. 2b) or fifth state (tree not shown). This is pre¬ 
sumably a result of putatively informative phylogenetic sig¬ 
nals contained within the larger ITS1 region of A. salina as 
opposed to all other species. A. urmiana and A. tibetiana 
are invariably grouped with the cluster of parthenogenetic 
populations, while A. sinica is more closely related to 


parthenogens from China (J1A7) than to all remaining 
Asian bisexuals and parthenogenetic strains. Within par¬ 
thenogens, three distinctive clusters are evident: the first 
one includes clones MAD4, TCL1 and TCL2 (closely relat¬ 
ed to A. urmiana ), the second one includes the other repre¬ 
sentative from Madagascar (MAD6) together with the two 
representatives from Namibia (NAM1 and NAM5) and 
one isolate from Greece (POL14), while the third cluster 
includes the remaining Greek isolates (MEM4 and 
MEM9) and parthenogens from Israel (EIL1 and EIL2). 

All previous groupings are faithfully retained irrespec¬ 
tive of the tree building method. Maximum likelihood 
and Bayesian inferred topologies are shown in Figs. 3 
and 4, respectively. Both methods yield a significantly bet¬ 
ter separation between A. franciscana from Great Salt Lake 
(GSL1 and GSL3) and the remaining sequences of this 
cluster compared with MP. All phylogenetic algorithms 
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Fig. 2. Majority rule (50%) maximum parsimony trees derived from the ITS1 sequence data showing phylogenetic relationships among Artemia species, 
(a) Gaps were treated as separate binary characters (see section 2.3) according to “simple indel coding” by Simmons and Ochoterena (2000). (b) Gaps were 
treated as missing data. Nodal support is indicated above branches. Bootstrap values (1000 replicates) of less than 60% for internal nodes are not shown. 
Thickened lines indicate parthenogenetic lineages. The topology is rooted with the outgroup Streptocephalus proboscideus. 
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Fig. 2 ( continued ) 


(MP, ML, and Bayesian) returned the same topology with 
A. persimilis as significantly differentiated from all the 
remaining species, A. franciscana invariably resolved as 
the sister group to Asian bisexuals and parthenogenetic 
Artemia and A. sinica as clearly divergent from the remain¬ 
ing Asian sexuals and the parthenogens (see Figs. 2a and b, 
3 and 4). 

The mtDNA phylogeny of a subset of taxa with corre¬ 
sponding ITS1 sequences is shown in Fig. 5. Although A. 
persimilis, A. franciscana, A. salina, and A. sinica are mark¬ 
edly differentiated, most parthenogens display identical 
mtDNA profiles with A. urmiana and A. tibetiana. Also, 
Greek parthenogenetic populations (POL and MEM) are 
grouped separately as they possess additional clones, two 
of which are private (D and G, see Table 3). Results of a 
Mantel test (Mantel, 1967) reveal a lack of association 


(r = 0.56, P > 0.05) between nuclear and mtDNA diver¬ 
gences for clones with complete ITS1 and 16S rRNA data. 
Thus, parthenogens with identical mtDNA clones (e.g. 
NAM, TCL, MAD, and EIL) are separated by 1TS1 diver¬ 
gences ranging from 0.5 to 3.9% (data not shown). The 
reverse is also true: the ITS1 distance between MEM and 
POL isolates is 0.6% while the haplotype divergences of 
their component clones range from 1.1 to 6.6% (Table 3). 

3.5. Age estimates 

Rate constancy for the entire ITS 1 data set was evaluat¬ 
ed following the method by Takezaki et al. (1995). The 
branch length test was applied to an initial neighbor joining 
topology, based on the HKY85 + G model of substitution. 
The assumption of constant rates of molecular evolution 
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0.1 

Fig. 3. Maximum likelihood phylogram (—In L = 6720.69) of Artemia species based on the HKY85 + G model of nucleotide substitution (—In 
L = 6752.37; A = 0.24, C = 0.20, G = 0.28, T = 0.28; number of substitution types Nst = 2; T s /T v = 0.99; y-shape parameter a = 1.94). Bootstrap values 
(100 replicates) of less than 60% for internal nodes are not shown. Boldface acronyms indicate parthenogenetic lineages. The topology is rooted with the 
outgroup Streptocephalus proboscideus. 


was rejected (U = 78.42, df =45, P < 0.05). Inspection of 
branch lengths resulted in exclusion of 10 sequences 
(WNA5, WNA4, EIL1, EIL2, SIN3, MAD6, ING3, 
TCL2, EIRM2, and TDP3) which invariably evolved signif¬ 
icantly different than the average rate. It has to be noted 
that, a number of these sequences had already been recog¬ 
nized as suspect outliers with the disparity index test (see 
section 3.1). The U statistic for the remaining sequences 
was 32.53 (P > 0.05, df = 35) and therefore, the hypothesis 
of clocklike evolution was justified. A linearized tree was 
then constructed for the pruned data set (Fig. 6). 

Estimates of times of divergence were based on a combi¬ 
nation of phylogenetic signal within the data matrix and 
major/regional biogeographic events. The former included 


the extreme differentiation and restricted presence of A. 
persimilis in Argentina and Chile, the clustering of A. salina 
into eastern and western Mediterranean groups, and the 
partitioning of the whole genus on the basis of ITS1 length. 
The latter were related to the breakup of Gondwana (see 
Nei, 1987; Mayr, 2001), to geological and hydrographic 
events in the Mediterranean basin during the Messinian 
salinity crisis (Penzo et al., 1998; Hrbek and Meyer, 2003; 
Huyse et al., 2004), and to the geological history of the 
Tibetan plateau. Hence, the calibration of A. persimilis 
divergence from the common ancestor was set to the time 
corresponding to the severance of South America from 
Gondwana at the middle Cretaceous (80-90 MYA). This 
yields a rate of sequence divergence of 0.2%/MYA within 
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Fig. 4. Bayesian inference topology of phylogenetic relationships in the genus Artemia. Sequence evolution is based on the HKY85 + G model. In two 
simultaneous runs, four Markov chains (three heated and one cold) were started from a random tree and run for 1,000,000 generations with sampling 
frequency every 100 generations. First 2000 trees were discarded as burn-in. Posterior probabilities higher than 0.90 are shown. Thickened lines indicate 
parthenogenetic lineages. The topology is rooted with the outgroup Streptocephalus proboscideus. 


each lineage. Artemia salina seems to have diverged ca. 
40 MYA, at about the time of emergence of the Mediterra¬ 
nean Sea, while its separation into eastern and western 
groups (ca. 7 MYA) roughly coincides with the beginning 
of the salinity crisis and subsequent fragmentation of the 
Mediterranean. Asian bisexuals and parthenogens last 
shared a common ancestor with A. franciscana 32 MYA, 
at the timeframe (early Tertiary) of a broad land connec¬ 
tion between North America and Europe (Mayr, 2001). 
Within the Asian clade of bisexuals, A. sinica branched 
off the remaining groups 8 MYA, within a period for which 
geological and thermochronological evidence (Clark et al., 
2005; Molnar, 2005) demarcate the elevation of the Tibetan 
plateau and subsequent separation of eastern and western 
landscapes. Considering our sampling effort, at least four 


independent origins of parthenogenesis are deduced, the 
most ancient of which dates at 3.5 MYA (A. sinica vs 
JIA7). The ages of the remaining parthenogenetic lineages 
could not be estimated since monophyletic arrangements 
were uncertain. 

4. Discussion 

The current study is the first robust, global phylogeny of 
the genus Artemia. Vicariance and pronounced ecological 
diversification seem to have shaped current biogeographic 
patterns to a great extent. The parthenogenetic gene pool 
is a heterogeneous amalgam of clones, including both nar¬ 
row endemics and widespread lineages of distinctive spatial 
and temporal origins and genetic diversity. Considering the 
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Fig. 5. NJ phylogram of Artemia populations assayed for 16S mtDNA 
RFLP patterns. The tree is based on nucleotide divergence estimates from 
restriction sites data. Asterisks indicate parthenogenetic populations. 


diversity of life cycles and modes of reproduction in zoo- 
planktonic invertebrates, our study provides additional 
insights into patterns of distribution in the continental 
aquatic realm and, in particular, for halophilic lineages 
that have made the transition to asexuality. It also corrob¬ 
orates the antiquity of anostracan members at the lower to 
middle Cretaceous (see Remigio and Hebert, 2000 and ref¬ 
erences therein; Daniels et al., 2004) and yields evidence 
that their historical biogeography may have been affected 
by common geological landmarks. 

The phylogenetic signal present in ITS 1 has proven suf¬ 
ficient in delineating Artemia species and, in some cases, 
detecting substantial differentiation at the intraspecific 
level. The observed interspecific divergence estimates (zero 
to 39.5%) are within the range reported for other continen¬ 
tal aquatic crustaceans of Gondwana origin. For example, 
in the fairy shrimp genus Streptocephalus, interspecific val¬ 
ues based on 5.8S-ITS1-18S sequence divergence vary 
from 0.90 to 90% (Daniels et al., 2004). However, lower 
ITS1 estimates, in the range of 15-21%, have been obtained 
for Darwinula (Schon et al., 1998). Regarding the intraspe¬ 
cific level, our estimate of 2.0% sequence diversity within 
parthenogenetic Artemia is striking when compared with 
Darwinula stevensoni. This ostracod is believed to be repro¬ 
ducing asexually for millions of years, however no ITS1 
variability was detected in a sampling range spanning from 
Finland to South Africa (Schon et al., 1998). A slow down 
in the rate of molecular evolution for ITS1 has been 
hypothesized to explain such extensive invariability. In 
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Fig. 6. Linearized tree of Artemia ITS1 sequences. Ten sequences evolving significantly different than the average rate at the 0.05 level were excluded 
(WNA5, WNA4, EIL1, EIL2, SIN3, MAD6, ING3, TCL2, URM2, and TDP3) and branch lengths were re-estimated following Takezaki et al. (1995). 
Triangles at the tips of the tree denote grouping of sequences for the purposes of illustration. The size of triangles is proportional to within-group 
divergence. Brackets indicate ITS1 length variation among clades. Values above and below the scale bar are MYA and sequence divergence, respectively. 
Parthenogens 1, POL14; NAM1; NAM5; Parthenogens 2, MEM4; MEM9; Parthenogens 3, MAD4; TCL1; and Parthenogens 4, JIA7. 
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our study, the obtained rate of 0.2%/MYA is generally 
higher than that of most Darwinulidae lineages with an 
exceptional fossil record (Schon et al., 2003) and may be 
a consequence of brine life as suggested by Hebert et al. 
(2002) for a number of halophilic lineages including Art- 
emia. Apart from its mode and rate of evolution, ITS1 
has been used for estimating phylogenies at different hierar¬ 
chical levels (Harris and Crandall, 2000; Weekers et al., 
2001; Chaw et al., 2005). Although ITS1 length variation 
has been reported in a number of organisms (Schlotterer 
et al., 1994; Von der Schulenburg et al., 2001), our results 
represent an unusual case of both increased size (compared 
with the ITS1 length of other anostracans, see Daniels 
et al., 2004) and clade-specific length variation within a sin¬ 
gle organismal lineage. Thus, apart from its phylogenetic 
utility, ITS1 may also serve as a crude marker for identifi¬ 
cation of unknown Artemia isolates. For example, isolates 
MAS1, MAS6, and MAS7 from Italy were confirmed as 
belonging to A. franciscana introduced samples (see Mura 
et al., 2004) while the species status (also A. franciscana) 
of the unknown ING3 isolate from Ingebright, Canada 
(see Table 1) was clarified. Furthermore, detailed investiga¬ 
tions of ITS 1 structure should provide useful insights into 
the origin, nature and possible association to high sequence 
change of ITS1 elements. 

The clade of A. persimilis diverged from the ancestral 
stock of Gondwana populations as a result of the tectonic 
separation of South America from Africa. The intraspecific 
component of diversity within A. persimilis is in agreement 
with previous analyses based on mtDNA markers (Gaj- 
ardo et al., 2004). The evolutionary route of A. salina seems 
to be closely related to the geologic history of the Mediter¬ 
ranean basin. The substantial differentiation into eastern 
and western A. salina groups, also strongly supported by 
AFLP analysis (Triantaphyllidis et al., 1997), is in accor¬ 
dance with the expansion of land and extensive desiccation 
of the Mediterranean Sea during the Messinian salinity cri¬ 
sis (see Penzo et al., 1998). The area of land exposure 
between southern Turkey and Egypt maps well with the 
eastern A. salina cluster, while that between the Gulf of 
Sidra and the Strait of Sicily coincides with the western 
A. salina cluster. Unlike previous studies (Abreu-Grobois, 
1987) indicating A. franciscana as the most variable and 
sub-structured species, our results point out that genetic 
diversity within A. salina may be at least 30% higher (see 
Table 2). Considering the island-like nature and abiotic 
variability of Artemia habitats, interpopulational differenc¬ 
es may be further enhanced and established rapidly, thus 
generating geographic population structure and leading 
to evolutionary groups with distinct ecological and life his¬ 
tory profiles (Kappas et al., 2004). 

Although the reconstructed phylogeny is concordant 
with well-established geological events, a degree of uncer¬ 
tainty exists over the topological placement of A. francis¬ 
cana , which is unconventional to current phylogenetic 
theory in Artemia. The affinity of A. franciscana with 
all Asian Artemia and parthenogens is geologically 


problematic, yet recurring. RAPD analysis (Sun et al., 
1999) has shown that the genetic distance between A. sinica 
and parthenogenetic Artemia is equal to that between A. 
franciscana and coastal Chinese parthenogens. Phylogenet¬ 
ic reconstruction based on COI and Cytb sequences (Perez 
et al., 1994) has indicated that A. franciscana shares a more 
recent relationship with parthenogenetic Artemia than with 
A. salina. More importantly, though, aspects of phenotypic 
evolution in the genus lend credence to the observed rela¬ 
tionship between A. franciscana and Asian-parthenogenetic 
lineages. Analysis of morphological characters (Trianta¬ 
phyllidis et al., 1997) has shown that a pair of characteristic 
spine-like outgrowths at the basal part of each penis is 
present in all bisexual Artemia except A. salina. This char¬ 
acter is also absent from male specimens of the extinct pop¬ 
ulation from Lymington, UK (collection of the Natural 
History Museum of London). Most remarkably, these 
structures are present in rare males ( Mura and Brecciaroli, 
2002) known to originate from asexual forms. The same 
groupings are also revealed by examination of the bulges 
at either side of the base of the ovisac (absent in A. salina ) 
(Amat, 1980) and the morphology of frontal knobs (sub- 
conical in A. salina as opposed to spherical in all other 
species) (Mura and Brecciaroli, 2004). These phenotypic 
patterns cannot be interpreted as primitive absences but 
rather as secondary losses or modifications in conjunction 
with a parallel ecological shift (low temperature tolerance, 
see Browne, 1992) in A. salina. The latter two characters, in 
particular, are relevant to reproductive morphology as they 
may determine the mechanical configuration of the couple 
during pairing. Our hypothesis and evolutionary timescales 
imply that the conflicting relationship between A. francis¬ 
cana and Asian-parthenogenetic lineages may be resolved 
by invoking phylogenetic explanations analogous to other 
zooplanktonic invertebrates (see Gomez et al., 2002; De 
Queiroz, 2005) and as yet unconsidered in Artemia. 
Although our results for A. franciscana may be regarded 
as provisional, they strongly suggest that current thinking 
on the primary separation of New and Old World lines 
may require a substantial re-evaluation. 

The present study highlights the benefits of Artemia as a 
model organism for investigating the interaction between 
sexuality and parthenogenesis. On the basis of nuclear 
DNA sequences and cytoplasmic markers, a genetically 
diverse assemblage of parthenogens is inferred, in close 
affinity with Asian bisexuals (A. sinica, A. urmiana, and 
A. tibetiana ). Although bootstrap values within the 
Asian-parthenogenetic clade are relatively wide for 1TS1 
(<60 to 100%), both nuclear and mtDNA phylogenies are 
consistent with a polyphyletic origin of parthenogens and 
suggest that hypotheses of a single origin (Beardmore 
and Abreu-Grobois, 1983) are incorrect. All Asian bisexu¬ 
als are genetically distinct, at least for nuclear DNA (see 
also Abatzopoulos et al., 2002b), and are implicated in 
transitions to asexuality. For example, the presence of a 
common haplotype (C, see Table 3) indicates that a num¬ 
ber of clones may have captured mtDNA from A. urmiana 
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or A. tibetiana. This implies that, for some genetic or eco¬ 
logical reasons, loss of sex is more frequent in the Old 
World. Due to the high nodal support from their closest 
sexual progenitor (A. sinica), only Chinese parthenogens 
(JIA7) could be confidently dated, bearing an age of 
3.5 MYA. Real antiquity can be debatable (Martens 
et al., 2003) and several authors (Law and Crespi, 2002 
and references therein) have proposed that generation time 
may be biologically more relevant than chronological time 
in assessing the evolutionary implications of asexual repro¬ 
duction. Taking into account that Artemia completes an 
average of 10 generations per year, a conservative estimate 
of 30,000,000 generations is obtained. Escaping mutational 
meltdown (Lynch et al., 1993) and extinction in such a gen¬ 
eration-span may be mediated by a combination of three 
mechanisms, all present in Artemia: large population sizes, 
protection provided by the dormant stage and polyploidy. 
It is therefore conceivable that apomictic clones of Artemia 
form monophyletic groups with their closest sexual ances¬ 
tors that are at least as old as the estimated divergence 
between A. sinica and the parthenogenetic population from 
China. A possibility exists, however, that this timeframe is 
not accurately estimated as sampling of Chinese partheno¬ 
gens was rather limited. Overestimation seems more unlike¬ 
ly though, since intraspecific diversity in A. sinica is known 
to be insignificant (Naihong et al., 2000) and at the order of 
magnitude determined here, unless of course a more recent 
sexual progenitor has gone extinct. 

A conclusion that can, with a certain degree of cautious¬ 
ness, be extended to all parthenogenetic representatives, 
regardless of their sources of origin and ages, is that their 
observed genetic diversity of 2.0% (see Table 2) is not a 
reflection of the diversity of their closest bisexual ancestors. 
At least for a common group of clones tested here, matrix 
comparison of nuclear and mtDNA divergences strength¬ 
ens the above reasoning. The lack of association between 
the diversities found within parthenogens and their nearest 
sexuals relatives, if verified, may provide clues to the mode 
of origin of parthenogenetic populations and designate 
intra- or interspecific hybridization as possible candidates 
(Schon et al., 2000; Delmotte et al., 2001, 2003). In the first 
case, the role of parthenogenetically produced rare males 
wifi be crucial and concrete data on their functionality 
and reproductive potential wifi be needed (Browne, 1992; 
Simon et al., 2003). On the other hand, interspecific hybrid¬ 
ization in the strict sense (crosses between bisexual species) 
would likely provide new insights into the distribution of 
parthenogenetic lineages, since their occurrence should, 
to a greater or a lesser extent, be associated with the distri¬ 
butions of their sexual parentals. In the long term, howev¬ 
er, physiological and ecological tolerances as well as genetic 
architecture (ploidy level, parthenogenetic mechanism of 
reproduction, i.e. apomixis or automixis) may also deter¬ 
mine niche breadth of asexuals (Zhang and Lefcort, 
1991). Although studies have reported marked ecological 
diversification among natural clones of Artemia (Browne 
and Hoopes, 1990; Baxevanis and Abatzopoulos, 2004), 


there is no evidence favouring either a general purpose or 
a specialist genotype. Our data also point to no specific 
direction, but a prediction made by Vrijenhoek (1998) for 
the possibility of intermediate clones occupying a narrow 
but universally available niche seems very promising in this 
case. This may also better explain competition between 
asexuals and their sexual relatives as well as the enigmatic 
absence, or at least unverified existence, of parthenogens 
from the New World. 
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